Derivation of the Maximum a Posterori Estimate for Discrete Time Descriptor Systems

In this report a derivation of the MAP state estimator objective function for general (possibly nonsquare) discrete time causal/non-causal descriptor systems is presented. The derivation made use of the Kronecker Canonical Transformation to extract the prior distribution on the descriptor state vector so that Maximum a Posteriori (MAP) point estimation can be used. The analysis indicates that the MAP estimate for index 1 causal descriptor systems does not require any model transformations and can be found recursively. Furthermore, if the descriptor system is of index 2 or higher and the noise free system is causal, then the MAP estimate can also be found recursively without model transformations provided that model causality is accounted for in designing the stochastic model. Index Terms Descriptor Systems, Maximum a Posteriori Estimate, Maximum Likelihood Estimate I. MAXIMUM A POSTERIORI ESTIMATION FOR DISCRETE TIME DESCRIPTOR SYSTEMS A. Introduction The first objective of this chapter is to find the maximum a posterior (MAP) estimate for the state vector sequence xk given a stochastic discrete time descriptor system model (SDTDS), noisy measurements and ∗Corresponding author. Email: [email protected] Monday 5th May, 2014 DRAFT 2 an informative prior for Ex0 as follows: Exk+1 =Axk +Buk + Fwk (1) yk =Hxk + vk (2) Ex0 ∼ N(r̄0, P0) (3) where N(r̄0, P0) denotes a normally distributed random variable with mean r̄0 and variance P0. We assume here that only the sequence uk is deterministic and all other sequences are random. The input disturbance sequence wk ∈ R and the measurement noise sequence vk ∈ R are i.i.d. normal random sequences; wk ∼ N(0, Iq) and vk ∼ N(0, R), where R ≻ 0. Furthermore, the random variables Ex0, wk, vk are assumed uncorrelated between each other. The matrices E,A ∈ Rneq×n, B ∈ Rneq×j , y ∈ R, H ∈ Rm×n and F ∈ Rneq×p. The maximum a posteriori estimate of xk is defined as the mode of the posterior distribution denoted by x̂ and given by: x̂ :=argmax x px|y(x|y) =argmax x py|x(y|x) px(x) py(y) = argmax x (log py|x(y|x) + log px(x)) (4) where, x̂ = {x̂ k } T k=0, x = {xk} T k=0, y = {yk} T k=0 and x̂ map k is the MAP estimate at time k. Linear descriptor systems define xk implicitly and hence the prior distribution px(x) can not be found directly from the stochastic descriptor system given in (1). A proceeding step is needed to convert the stochastic descriptor system to a format that reveals the prior distribution on xk. On the other hand, the constrained maximum likelihood estimate x̂ k is found by treating the state sequence xk as a parameter and the estimates are obtained by maximizing the likelihood function: x̂ k := argmax x L(x|y) = argmax x py|x(y|x) subject to x ∈ C (5) where the constraint x ∈ C forms the prior information about the parameter x. In [1], the input sequence Buk and the prior for Ex0 were reformulated as noisy measurements: Buk =Exk+1 −Axk − Fwk r̄0 =Ex0 + e where e is a Gaussian zero mean random vector with variance P0 and independent of wk, vk, while xk was viewed as parameters. The objective was to construct recursively the filtered or predicted estimate Monday 5th May, 2014 DRAFT 3 given by the conditional mean. It may be argued, however, that this paradigm shift in viewing Buk as a measurement is inconsistent with the reality that Buk is a user defined input that is not random. Also, the study in [2] presented an algorithm for transforming non-causal stochastic descriptor systems into causal systems but did not analyse how to avoid stochastic non-causality which is more meaningful for state estimation problems in practice. In [3] and [4], matrix and state variable transformations were used to recast state estimation problems for square causal/non-causal descriptor systems into conventional state space estimation problems. However, the method results in estimating transformed state variables instead of the original model variables and hence adds a requirement for an inverse transformation at every iteration for finding the estimates. Moreover, the method was not generalized to non-square descriptor systems. In this chapter, it is shown that model transformations are not necessary if the system is causal and the algebraic equations are modelled properly to avoid stochastic non-causality. The analysis is based on examining the various subsystems that descriptor systems can represent using Kronecker canonical transformation. This canonical form is suitable for extracting the prior on xk for the most general case of the system dynamics (1) (i.e. causal or non-causal, square or non-square), and is also capable of revealing the necessary assumptions and restrictions needed on the stochastic model and noisy measurements that define a well posed estimation problem. B. The Real Kronecker Canonical Form of a Matrix Pencil λE −A The Kronecker canonical form transformation (KCF) for singular matrix pencils λE−A was developed by the German mathematician Leopold Kronecker in 1890. This is also often called the generalized Schur decomposition of an arbitrary matrix pencil λE−A and is a generalization of the Jordan canonical form for a square matrix. Definition I.1. [5] The matrix pencil λE−A is said to be singular if neq 6= n or det(λE−A) = 0 ∀λ ∈ C. Otherwise, if neq = n and there exist a λ ∈ C such that det(λE − A) 6= 0 then the matrix pencil is called regular. Definition I.2. [5] The matrix pencil λẼ− Ã is said to be strictly equivalent to the matrix pencil λE−A for all λ ∈ C if there exist constant non-singular matrices P ∈ Cneq×neq and Q ∈ Cn×n independent of λ such that: Ẽ = PEQ, Ã = PAQ (6) Monday 5th May, 2014 DRAFT 4 Definition I.3. A nonzero vector x ∈ C is a generalized eigenvector of the pair (E,A) if there exists a scalar λ ∈ C, called a generalized eigenvalue such that: (λE −A)x = 0 Theorem I.4. [5], [6] Let E,A ∈ Rneq×n. Then there exists non-singular matrices P ∈ Rneq×neq and Q ∈ Rn×n for all λ ∈ R such that: P (λE −A)Q = λẼ − Ã = diag ( Uǫ0 , · · · , Uǫp , Jρ1 , · · · , Jρr , Nσ1 , · · · , Nσo , Oη0 , · · · , Oηq ) (7) where the matrix blocks are defined as follows: 1) The block Uǫ0 correspond to the existence of scalar dependencies between the columns of λE − A and is a zero matrix of size neq× ǫ0. Blocks of the type Uǫi for i = 1, · · · p are the bidiagonal pencil blocks of size ǫi × (ǫi + 1) and have the form: Uǫi = λEUǫi −AUǫi = λ 